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' Abstract. Thermodynamic properties of charge-stabihsed colloidal suspensions and 

04 ■ polyelectrolyte solutions are commonly modeled by implementing the mean-field 

' Poisson-Boltzmann (PB) theory within a cell model. This approach models a bulk 

f"]^ ' system by a single macroion, together with counterions and salt ions, confined to 

O . a symmetrically shaped, electroneutral cell. While easing numerical solution of the 

^ • I nonlinear PB equation, the cell model neglects microion-induced interactions and 

' correlations between macroions, precluding modeling of macroion ordering phenomena. 

^ , An alternative approach, which avoids the artificial constraints of cell geometry, 

^ ' exploits the mapping of a macroion-microion mixture onto a one-component model of 

pH . pseudo-macroions governed by effective interparticle interactions. In practice, effective- 

O I interaction models are usually based on linear-screening approximations, which can 

I accurately describe strong nonlinear screening only by incorporating an effective 

(renormalized) macroion charge. Combining charge renormalization and linearized 
I PB theories, in both the cell model and an effective-interaction (cell-free) model, 

«) ■ we compute osmotic pressures of highly charged colloids and monovalent microions, 

7-H , in Donnan equilibrium with a salt reservoir, over a range of concentrations. By 

comparing predictions with primitive model simulation data for salt-free suspensions, 
and with predictions of nonlinear PB theory for salty suspensions, we chart the 
"nI I limits of both the cell model and linear-screening approximations in modeling bulk 

■ thermodynamic properties. Up to moderately strong electrostatic couplings, the 

cell model proves accurate in predicting osmotic pressures of deionized (counterion- 
dominated) suspensions. With increasing salt concentration, however, the relative 
contribution of macroion interactions to the osmotic pressure grows, leading predictions 
of the cell and effective-interaction models to deviate. No evidence is found for a 



o 
o 



X 

' liquid-vapour phase instability driven by monovalent microions. These results may 

guide applications of PB theory to colloidal suspensions and other soft materials. 
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1. Introduction 

When dispersed in water, or other polar solvents, colloidal particles can acquire an 
electric charge through surface dissociation of counterions Electrostatic repulsions 
between charged macroions, screened by surrounding microions (counterions, salt ions), 
act to stabilise colloidal suspensions against aggregation induced by attractive van der 
Waals forces, as first explained by the classic theory of Derjaguin, Landau, Verwey, 
and Overbeek (DLVO) [2[|3]. Equilibrium and dynamical properties of many soft and 
biologically prevalent materials, including charge-stabilised colloids, polyelectrolytes, 
and ionic surfactants, are largely determined by microion distributions and induced 
electrostatic interactions between macroions. 

Molecular simulations provide essentially exact results for model colloidal 
suspensions and polyelectrolyte solutions, but can access only limited parameter ranges 
of these complex systems. For salt-free suspensions, recent Monte Carlo studies of the 
primitive model have probed osmotic pressure, thermodynamic stability, and macroion 
structure [U O |6l [TJ El [9] and characterized effective macroion interactions [TOl [TTl 
[T2I [T3I [T4] and bulk phase behaviour, including liquid-vapour separation [151 IlHl [12] 
and melting [18]. At significant salt concentrations, proliferation of microions renders 
simulations of the primitive model unwieldy, and data for bulk properties of salty 
suspensions are scarce. 

Materials with strong electrostatic interparticle interactions ^2] are commonly 
modeled via Poisson-Boltzmann (PB) theory [201 El], which combines the exact Poisson 
equation for the electrostatic potential with a mean-field (Boltzmann) approximation 
for the microion densities. Neglect of correlations between microions is strictly 
justifiable only for weak electrostatic couplings, characteristic of monovalent microions 
in water. Nevertheless, the mean-field PB theory often makes at least qualitatively 
correct predictions for thermodynamic properties and provides a valuable reference for 
simulations and more sophisticated theories that incorporate correlations. 

Central to Poisson-Boltzmann theory is the PB equation governing the electrostatic 
potential. For symmetric (isotropic) boundary conditions, this nonlinear differential 
equation is easily solved — analytically in planar and cylindrical geometries, and 
numerically in spherical coordinates. For this reason, PB theory is most often 
implemented within a cell model [2T| [22] , which reduces a bulk suspension to a single 
macroion — idealised as a charged sphere, rod, or plate — centred in a cell of the same 
shape, together with salt ions and neutralizing counterions. While a symmetric cell may 
reasonably approximate the Wigner-Seitz cell of a periodic crystal, the computational 
advantage gained by imposing symmetry is offset by a loss of accuracy in neglecting 
microion-induced interactions and correlations between macroions. Although efficient 
for modeling nonlinear microion screening, the cell model is not designed to describe 
macroion structure and ordering. 

An alternative theoretical approach, which strives to incorporate macroion 
structure, maps the multi-component ion mixture onto a one-component model by 
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averaging over microion degrees of freedom in the partition function [23] • The 
resulting pseudo-macroions are governed by an effective Hamiltonian involving effective 
electrostatic interactions. Although the one-component model is free of the constraints 
of cell geometry, practical implementation usually requires invoking linear-screening 
approximations, derived by linearizing the PB equation. Moreover, in applying 
linearized theories to strongly-coupled suspensions of highly charged macroions, it proves 
essential to consistently incorporate an effective (renormalized) macroion charge to 
properly account for nonlinear microion screening [211 ESI EHJ [2Zj • 

Condensation of counterions onto highly charged polyelectrolyte chains is well 
established [2HI EH]. By comparison, the association of counterions with spherical 
macroions is a rather more subtle and less understood phenomenon. Several theories 
have been proposed to describe counterion accumulation around highly charged 
macroions and the resulting reduction of the effective macroion valence. In a pioneering 
study, Alexander et al. [30] defined an effective macroion valence for the Debye-Hiickel 
(linear-screening) theory by matching solutions of the nonlinear and linearized PB 
equations at the boundary of a spherical cell. In a series of studies. Levin and 
coworkers applied liquid-state theory [3ll [321 EH] and derived an effective valence from 
a jellium approximation [SlJ |35l [361 EZ] fo^' the total free energy. Other thermodynamic 
approaches are based on a thermal criterion for either the electrostatic potential around 
a macroion [38l|39lll0] or the effective pair potential [H]. Recently, an effective valence 
has been incorporated into the one-component model with linear-screening effective 
interactions [121 Sal HI US] . 

Effective macroion charges have been measured in experiments, e.g., by torsional 
resonance spectroscopy [46], light or small angle neutron scattering [471 SH], and 
conductivity [19] measurements. Renormalized charges also have been extracted from 
Monte Carlo (MC) simulations according to various criteria: equating the Debye-Hiickel 
pressure [TD] or chemical potential [50] to the corresponding MC quantity; applying 
the prescription of Alexander et al. [30] to the spherical cell model [6]; fitting an 
effective Yukawa pair potential by inverting the macroion-macroion pair distribution 
function [13]; from the inflection point of the curve of accumulated running charge vs. 
inverse separation from the macroion centre [HI |51]; and from a dynamical criterion 
relating the counterion kinetic and potential energies [52] . 

The purpose of this paper is twofold: first, to generalize the charge renormalization 
theory proposed in ref. [S] from salt-free suspensions to suspensions dialyzed against an 
electrolyte reservoir; second, to apply the generalized theory to salty suspensions to test 
the accuracy of the PB cell model in predicting osmotic pressure over a wide range of 
salt concentrations. By comparing results with available simulation data for deionized 
suspensions and with predictions of nonlinear PB theory for salty suspensions, we 
explore the limits of both the cell model and linearization approximations in predicting 
bulk thermodynamic properties. Up to moderately strong electrostatic couplings, the 
cell model proves accurate in predicting osmotic pressures of counterion-dominated 
suspensions (low salt concentrations). With increasing salt concentration, as macroion 
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interactions gain in relative importance, predictions of the one- component (effective- 
interaction) model increasingly deviate from those of the cell model. No evidence is 
found of a liquid-vapour spinodal instability for a physically rational definition of the 
renormalization threshold. 

The paper is organized as follows. Section |2] first reviews the relevant models: the 
primitive model of charged colloids, the charge renormalization model for the effective 
macroion valence, and the one-component model governed by an effective Hamiltonian. 
Working from the density-functional formulation of Poisson-Boltzmann theory. Sect. [3] 
then derives the PB equation and outlines the two main implementations: the cell 
model and the effective- interaction model. In Sect. HI the two implementations of PB 
theory are linearized about the mean potential or, equivalently, the mean microion 
densities. Section [5] applies the nonlinear and linearized PB theory to calculate the 
osmotic pressure and compares predictions with available data from primitive model 
simulations of salt-free suspensions. Finally, Sect. [6] summarizes and concludes. 

2. Models 

2.1. Primitive Model 

The primitive model of charged colloids and polyelectrolytes ideahses the solvent as 
a homogeneous dielectric continuum of relative permittivity e. Dispersed throughout 
the solvent are macroions and microions, here modeled, respectively, as charged hard 
spheres of radius a and valence Z (charge —Ze) and monovalent point ions. In a closed 
suspension, all particles are confined to the same volume V . In Donnan equilibrium, 
only the macroions are confined, while the microions can exchange (e.g., across a 
semi-permeable membrane) with a salt reservoir, here assumed to be a 1:1 electrolyte. 
Modeling the reservoir as an ideal gas of ions of pair number density rio, the suspension 
has fixed salt chemical potential /i^ = 2/xo, where /^o = /c_Bnn(noA^) is the chemical 
potential of each microion species at absolute temperature T, the thermal wavelength 
A defining the arbitrary zero of the chemical potential. 

A bulk suspension of Nm macroions, counterions, and A^^ dissociated pairs of 
oppositely charged salt ions contains A^+ = N^. + Ng positive and A^_ = Ng negative 
microions, for a total of = Nc + 2Ns microions. Global electroneutrality constrains 
the macroion and counterion numbers by the condition ZN^ = N^. Accounting for 
excluded-volume and electrostatic (Coulomb) pairwise interparticle interactions, the 
primitive model Hamiltonian can be expressed as 




(1) 



which separates naturally into three terms: a macroion Hamiltonian, 

Z^e^ ^ 1 





(2) 
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where i^hs is the macroion hard-sphere Hamiltonian (including kinetic energy) and 
is the centre-centre separation between ions; a microion Hamiltonian, 

2 

H, = K, + '-y:"-f^^ (3) 

with microion valences Zi = ±1 and kinetic energy K^] and a macroion-microion 
interaction energy, 

^^. = ^EE^- (4) 

2.2. Charge Renormalization Model for Effective Charge 

The parameter Z in the primitive model corresponds to the maximum (bare) charge 
of the macroions, measureable by titration. In practice, this simple correspondence 
may be complicated by chemical reactions at the macroion surface, dependent on pH 
and salinity, which can affect ion dissociation |171 |53l |5l]. Even in the absence of 
charge regulation, however, the effective macroion charge — probed by electrophoresis, 
rheology, and scattering experiments — may lie considerably below the bare charge. 

The distinction between bare and effective charges arises from the strong association 
between counterions and highly charged macroions. Planar and rodlike macroions 
generate electrostatic potentials of sufficiently long range to condense counterions above 
a threshold surface charge density [2^. In contrast, the potential outside of a spherical 
macroion decays too rapidly to overcome entropy and irreversibly bind counterions. 
Nevertheless, counterions that venture sufficiently close to the macroion surface to 
become thermally (transiently) bound can be considered to renormalize the bare charge 
to a lower effective charge [30] . The remaining, weakly associated, counterions contribute 
to screening of the effective charge and are amenable to linear-screening approximations. 

Akin to Oosawa's two-phase theory of polyelectrolyte solutions [29], free and bound 
counterions part ways at a distance from the macroion at which the deviation of the 
electrostatic energy of counterion-macroion attraction from the average potential energy 
rivals the typical thermal energy fc^T per counterion. Defining %l){r) as the reduced 
electrostatic potential energy (in ksT units) at distance r from the centre of a macroion, 
the thickness of the spherical association shell 5 (Fig. [1]) is defined by the condition 

|^(a + 5)-^|=C, (5) 

where ip is the spatially-averaged potential and C is a constant of 0{1) yet to be 
specified. At distances r > a + 5, the deviation is presumed to be sufficiently small 
that a linear-screening approximation is valid. Assuming colons to be expelled from the 
association shell, the "dressed" macroion [55], consisting of a bare macroion and its shell 
of bound counterions, carries an effective valence Z < Z. Statistical fluctuations in Z are 
neglected in this simple model. Other workers have applied a thermal criterion similar 
to Eq. ([5]) to the electrostatic potential [381 ESI HQ] or the effective pair potential [H] . 
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2.3. One- Component Model and Effective Hamiltonian 

Simulations of the primitive model [3 E] have characterized bulk properties of salt-free 
(counterion-dominated) suspensions with relatively low charge asymmetry [Z < 100). 
Significant concentrations of salt pose, however, computational challenges for large- 
scale simulations of bulk suspensions. Consider, for example, that at 10 mM salt 
concentration a suspension of 100 macroions of radius a = 10 nm at 1% volume fraction 
contains (9(10^) particles, all interacting via long-range Coulomb forces, and that the 
particle number scales as a?. Therefore, salt-dominated suspensions usually are modeled 
by first mapping the mult i- component mixture onto a one-component model (OCM). 
In Donnan equilibrium, a suspension is governed by a semigrand partition function, 

Z = {{expi-(3H))^)^ , (6) 

where ( )^ denotes a grand canonical trace over microion coordinates, ( )^ a canonical 
trace over macroion coordinates, and (3 = l/{kBT). The one-component mapping 
expresses the partition function in the form 

Z = (exp(-/3i/eff))^ , (7) 

where the effective Hamiltonian for a one-component system of pseudo-macroions 

H,s = Hm- kBT\ji{exp[-f3{H^ + H^^)])^ = (8) 

can be interpreted equivalently as the grand potential fi^ of the microions in the 
"external" potential of the macroions. For a closed suspension, with fixed microion 
numbers, ( )^ becomes a canonical trace over microion coordinates and is replaced 
by the Helmholtz free energy of the microions. 

Practical applications of the OCM require approximating Q^. For this purpose, 
Poisson-Boltzmann theory is a powerful approach. Below we briefly review PB theory 
and two common implementations: the cell model and the effective- interaction approach. 



3. Poisson-Boltzmann Theory 

3.1. Density- Functional Formulation 

Poisson-Boltzmann theory is most elegantly formulated within the framework of classical 
density-functional theory of nonuniform fluids [56l [57] . Corresponding to the primitive 
model Hamiltonian [Eqs. (l2])-(j4])], there exists a Helmholtz free energy functional 
F[nm{r),n±{r)], which (for a given external potential) is a unique functional of the 
macroion and microion number density profiles, nm{r) and n±{r), varying with spatial 
position r. This free energy functional separates, according to F = Fj^ + F^^ + F^^t, into 
a (purely entropic) ideal-gas free energy functional of all ions, 

Fid = kBT J dr n.{r){\n[n,{r)A'] - 1} , (9) 

an excess free energy functional, Fgx = Fhs+Fgi, due to hard-sphere (hs) and electrostatic 
(el) interparticle interactions, and a contribution Fext due to an external potential. 
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Neglecting interparticle correlations (mean- field approximation), the electrostatic part 
of the excess free energy functional may be approximated as 

i^ei = ^/^drp(r)$(r) , (10) 

where p(r) = e[n^{j) — n__(r) — n/(r)] is the total charge density, including the negative 
charge fixed on the macroion surfaces of number density ra/(r), and 

<l>(r) = / dr'-^^^ (11) 



e|r — r'l 

is the total electrostatic potential at position r due to all ions. The potential and the 
total ion density are related via the exact Poisson equation 

V2$(r) = -^p(r) , (12) 

which may be expressed in the form 

VV(r) = -47rABK(r) - n_(r)] ; V^|surface = ZXe/a^ , (13) 

where if) = /3e$ is the reduced electrostatic potential, = /3e^/e defines the Bjerrum 
length, and the macroion charges are absorbed into a boundary condition at the 
macroion surfaces. The microion densities implicitly vanish inside the macroion cores. 

For a given external potential, the equilibrium densities of all ions minimize the total 
grand potential functional of the system. Alternatively, fixing the macroion coordinates 
and regarding their charges as the source of an external potential, the equilibrium 
microion densities alone minimize the microion grand potential functional 

n^[n±{v)] = F^[n±{v)] - j^drn+{Y) - j^dvn^{v) , (14) 

defined as a Legendre transform of the microion free energy functional 

F^[n±{r)] = knT f dr J] n,(r){lnK(r)A V 1} + ^ / dr p(r)<l>(r) , (15) 

the microion (electro) chemical potentials fi± being identified as the Legendre variables. 
Note that fi^ depends parametrically on the macroion coordinates and that includes 
macroion-macroion Coulomb interactions for electroneutrality. Under the assumption 
that either the electrostatic potential or the electric field vanishes everywhere on the 
boundary of the volume V, the microion free energy functional also may be expressed 
in the form 

/3F^[n±(r)] = / dr ^ n,{r){ln[n,{r)A^] - 1} + -i— / dr |V^P . (16) 
Jv ottab Jv 

Minimizing fi^[n±(r)] with respect to n±{r) now yields the mean-field Boltzmann 
approximation for the equilibrium microion densities 

n±(r) = exp [=!='?/' (r)] (fixed macroions) , (17) 

the reference densities, ni°^ = A-3exp(/3/i±), being the microion densities at the 
reference potential ip = 0. Evaluating the microion grand potential functional [Eq. ( IT^ ] 
at the equilibrium microion density profiles [Eq. f fT7|) ] yields the microion grand potential 

/3Qf, = -f dr[n+(r)+n_(r)]-^ /" dr [n+(r)-n_(r)+n/(r)]^(r) . (18) 
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For a closed suspension (fixed particle numbers), the chemical potentials of the two 
microion species differ because of asymmetric interactions with the macroions: /i+ 7^ 
Correspondingly, the reference densities also differ: n+ ^ 7^ ul^ . In Donnan equilibrium, 
however, exchange of microions with a salt reservoir shifts the intrinsic microion chemical 
potentials, /x™ = [5F^/5n-|-(r)]eq, by the Donnan potential ipD'- 

/3/i± = ±i^D = /3/io = In(noA') . (19) 

The total chemical potentials, and so too the reference densities, of the two microion 
species are thus equalized. The equilibrium microion density profiles are then given by 

n±(r) = raoexp[=FV'(r)] . (20) 

The Donnan potential is interpreted physically as the change in electrostatic potential 
across the reservoir-suspension interface, and mathematically as a Lagrange multiplier 
for the constraint of global electroneutrality. 

Combining the Poisson equation for the potential [Eq. ( fT3ll ] with the Boltzmann 
approximation for the microion densities [Eq. (fT7|) or (!20l) ]. the Poisson-Boltzmann 
equation takes the form 

- I (^+^ exp[-^(r)] - exp[^(r)]) (closed) ^^^^ 

\ nlsmh'il){Y) ; Wlsurface = ^Afi/a^ , (Donnan) 

where = ^/^t^XbUq is the screening constant in the reservoir. Beyond the boundary 
condition at the macroion surfaces, the boundary-value problem is fully specified only 
by imposing another condition at the outer boundary of the system, which depends on 
the practical implementation of the theory. 



3.2. Cell-Model Implementation 

The anisotropic boundary conditions on the nonlinear PB equation imposed by 
an arbitrary configuration of macroions render a general solution of Eq. fl2Tl) 
computationally daunting. In recent years, powerful ah initio methods have been 
developed for combining PB theory of microion density profiles with molecular 
dynamics [56l [571 [58l [59] or Brownian dynamics [6OI [6T| simulation to evolve macroion 
coordinates according to derived forces. Despite such advances, however, most 
applications of PB theory have been implemented within a cell model to facilitate 
numerical solution. 

In a seemingly bold reduction, the cell model represents a bulk suspension by 
a single macroion, neutralizing counterions, and salt ions confined to a cell of the 
same shape as the macroion. Microion-induced interactions between macroions are 
simply ignored. For spherical colloids, the natural choice is a spherical cell centred 
on the macroion [62]. Gauss's law then dictates that the electric field must vanish 
everywhere on the boundary of the electroneutral cell. With the potential and microion 
densities depending on only the radial coordinate r, the PB equation reduces to an 
ordinary differential equation for %lj{r) with boundary conditions ilj'{a) = ZXe/o^ and 
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ip'{R) = 0, where the cell radius R is commensurate with the average macroion density: 
rim = N„JV = 3/{4ttR^). For a closed suspension, the arbitrary location of the reference 
point of the electrostatic potential (where = 0) is usually chosen as the cell boundary: 
ip{R) = 0. In Donnan equilibrium, the potential is conventionally chosen to vanish in 
the reservoir, in which case the boundary value of the electrostatic potential is identified 
as the Donnan potential: ip{R) = ipo 7^ 0. 

An appealing feature of the cell model is the simple analytic relation between the 
bulk pressure p and the microion densities at the cell boundary: 

= n+{R) + n_{R) . (22) 

Although first derived within the mean- field PB framework [22j, this pressure relation 
proves to be exact within the cell model [63], i.e., valid also for correlated microions. 
In Donnan equilibrium with an ideal-gas reservoir, the osmotic pressure 11, i.e., the 
difference in pressure between suspension and reservoir, is then given by 

/3n = n+{R) + n^{R) - 2no . (23) 



Within PB theory, the osmotic pressure is strictly positive [26], which follows simply 
from Eqs. ( l20i) and ( l23l) and the inequality, coshx > 1. 

The cell model implements PB theory by approximating the microion grand 
potential [Eq. ([8])] in the one-component mapping of the primitive model (Sec. 12. 3p . 
Reducing a bulk suspension to a single macroion in a symmetrically shaped cell with 
isotropic boundary conditions greatly facilitates solution of the nonlinear PB equation. 
The microion grand potential implicitly depends on the average macroion density 
through the density-dependent cell radius. Moreover, bare Coulomb interactions among 
all macroions are included in an average sense via the electroneutrality constraint. 
The cell-model implementation excellently approximates the single-macroion (one-body) 
contribution to the free energy and pressure. The trade-off for fully incorporating 
nonlinear microion screening, however, is complete neglect of effective interactions and 
correlations among macroions induced by microions. Accurately accounting for multi- 
macroion contributions to thermodynamic properties requires an effective-interaction 
implementation of PB theory. 



3.3. Effective-Interaction Implementation 

An alternative implementation of PB theory, also based on the one-component mapping, 
focuses on effective interactions derived from perturbative expansion of the microion 
grand potential about a uniform reference system, namely, a plasma of microions 
unperturbed by the macroions. By incorporating macroion interactions, this approach 
can model both thermodynamic and structural properties of colloidal suspensions. 
Several liquid-state theoretical frameworks have been developed. Density-functional 
theories [HH [65l [66l [67] expand the ideal-gas free energy functional in a Taylor series in 
powers of deviations of the microion density profiles from their mean values Uj-: 

/3Fid[n±(r)] = ^ (^N^nin^A^) - 1] + ^ J dr ^(r) - n,]' + • ■ ■) . (24) 
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Response theories [HI ESI ESI EH ESI ES] expand the microion density profiles about a 
reference plasma in powers of the macroion potential, 0m±(r) = / dr'fm±(|r — r'|)nm(r'): 



where Xij = '^^i(r)/50mj(r') are linear response functions of the reference plasma and 
higher-order terms involve nonlinear response functions. Distribution function theories 
are based on density-functional expansions of the microion correlation functions about 
a reference plasma [681 EH Effl ED E21 ES] or on integral-equation closures [801 El 182| 

EallElllESlES]. 

In all of these statistical mechanical schemes, insertion of the perturbative expansion 
into Eq. (jS]) recasts the effective Hamiltonian (microion grand potential) in the form of 
a sum of effective interactions: 



where is a one-body volume energy - independent of the macroion coordinates but 
dependent on the average macroion density - fefil'") is an effective electrostatic pair 
potential summed over macroion pairs with centre-centre separation r, and higher-order 
terms involve sums over effective many-body interactions. The volume energy here 
corresponds to the microion grand potential for a single macroion. 

In practice, mean-field (Debye-Hiickel-like) approximations for the electrostatic 
free energy, response functions, or correlation functions, make the various effective- 
interaction approaches essentially equivalent. Thermodynamic and structural properties 
of a bulk system can be calculated by inputting the effective interactions into statistical 
mechanical theories or simulations of the OCM. 

4. Linearized Poisson-Boltzmann Theory 
4-1. Perturbation Expansion and Truncation 

Practical applications of PB theory often invoke a perturbation approximation, whereby 
the microion densities on the right side of the PB equation [Eq. ( 12T|) ] are expanded in 
powers of the deviation of the potential from a reference potential. Although not unique, 
the reference potential is commonly chosen as the reservoir potential [ip = 0) for a 
suspension in Donnan equilibrium. A more appropriate choice, however, is the mean 
potential in the suspension ip, this being also the only consistent choice for a closed 
suspension [211 |25l |26l |271 El] . The distinction between these two choices is especially 
relevant for deionized suspensions of highly charged colloids. 

A linear-screening approximation truncates the expansions of the microion densities 
at first (linear) order in the potential. Correspondingly, the microion grand potential 
functional is expanded to quadratic order in the deviations of the microion densities 



ni{r) = ni+Y, / dr'xij(|r - r'|)0^j(r') H , 



(25) 




(26) 
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from their means 

/3fiJ"[n±(r)] = Y.\N. 

i=± ^ 

^ y dr [n+(r) - n_(r) - nf{v)]ip{r) 



ln(^ 

no 



1 

2ni 



dr hiir) 



+ 



(27) 



where n± = N±/V' here represent the average microion densities in the free volume V 
(unoccupied by the macroion hard cores). Minimizing the approximate functional with 
respect to n±(r) now yields the hnearized microion density profiles [of. Eq. fl20l) ]: 



ln(^U "^^^^-"^ ±^(r) 
, no / n± 







f28) 



and thus 



n±(r) =n±[lT^(r)±V^] , (29) 

where n± = riQexp^^ip) for consistency. Substituting Eq. f l29|) back into Eq. fl27j) yields 
the linearized microion grand potential per macroion: 



/3w 



lin 



E 



ln( — 

no 



(30) 



where x± = N±/Nm- The first term on the right side originates from the microion 
entropy; the second term combines the macroion self energy, Z^AB/(2a), and the 
macroion-counterion interaction energy; and the final term accounts for the average 
potential energy of the microions relative to the reservoir. Practical implementations 
aim to calculate the potential. 



4-2. Cell- Model Implementation of Linearized Theory 

Substituting the hnearized microion densities [Eq. ( 129|) ] into [Eq. ( 121]) ]. the linearized 
PB equation in spherical coordinates takes the form 

2 — 

^"{r) + -i)'{r) = K^[i){r)-^]-A-n\B{n+-n^) , a<r<R, (31) 



where k = y 47rAB(n+ + n_) is the screening constant in the suspension. Note that 
K depends implicitly on the average macroion density via the global constraint of 
electroneutrality and that cosh^ in Donnan equilibrium, but not in closed 

equilibrium. With cell boundary conditions, ip'la) = ZXb/o^ and ip'iR) = 0, the 
solution of Eq. f l3T]) is 



C 6 

ij{r) = Ai— + + As , a<r<R, 



where the coefficients are given by 



Ai = ZXb 



[na - l)e'"' - {na + 1) 



kR 



kR+1 



^K(2R-a) 



A2 = 



kR + 1 



(32) 

(33) 
(34) 
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A3 = ^ + . 35 

For a closed suspension, the choice ip{R) = gives a mean potential 

« = -^^e--^. (36) 

Kit +1 n+ + ?7,_ 

For a Donnan suspension (with zero reservoir potential), averaging Eq. (1321) over the 
free volume of the spherical cell yields rather 

^ = ^3 3 / ^"^ r'^^ir) = - sinh"^ -— ^ . (37) 

The pressure relation [Eq. ( |22l) ]. although exact in the cell- model implementation of 
nonlinear PB theory, does not survive linearization. In linearized PB theory, therefore, 
the pressure must be computed as a thermodynamic derivative of the microion grand 
potential [Eqs. ([27D and (|30l)]: 

\ / T,Nm,no \ / T,no 

4-3. Effective-Interaction Implementation of Linearized Theory 

In contrast to the cell model, effective-interaction implementations of PB theory include 
microion-induced interactions between macroions in the suspension. Truncating the 
ideal-gas free energy functional expansion [Eq. ( 12^ ] at quadratic order, or the microion 
density profile expansions [Eq. (125!) ] at linear order, amounts to neglecting three-body 
and all higher-order effective interactions. At the mean-field level, such linear-screening 
approximations predict an effective electrostatic pair potential of Yukawa form, 

Mr) = Z^Xb — , r > 2a , (39) 

\1 + Ka/ r 

and a one-body volume energy [E^ in Eq. ( I26l) ] identical in form to the linearized 
microion grand potential [Eq. fl30|) ]. The electrostatic potential differs, however, from 
the cell model prediction. In contrast to the solution of the linearized PB equation 
in the spherical cell model [Eqs. (I32l) -( l35|) ]. hnear-response theory [751 EHl EH] and 
related effective-interaction theories [Ml ESI [66] predict 

^(r) = -Z\b — , r > a , (40) 

1 + r 

with mean value 

for macroion volume fraction 77 = (47r/3)a'^(A^m/^) = 1 — V' /V . This potential is 
identical to the solution of the linearized PB equation for a single macroion in a bulk 
suspension under open boundary conditions: ?/'(r) and %l)'{r) — )■ as r — )■ 00. 
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Taking care to consistently incorporate density dependence, the effective electro- 
static interactions (including the volume energy) may be input into simulations [87j or 
statistical mechanical theories for the total semigrand potential 

n = -ksT \nZ = E^ + Fm . (42) 

In linear-response theory, the macroion free energy 

Fm = -ksTln (exp (^-f3H^, - ^ ^ f3v,s{r,,)^ ^ , (43) 

which includes contributions from the macroion entropy and interactions, can be 
reasonably approximated by a simple variational approach. First-order thermodynamic 
perturbation theory with a hard-sphere reference system [88] predicts a macroion free 
energy per macroion 

= min \fhs{d) + 27™^ / drr'^gh^ir; d)ves{r)\ , (44) 

(d) I Jd J 

where /hs and gi^s are the free energy per macroion and radial distribution function, 
respectively, of the reference system. Minimization with respect to the effective hard- 
sphere diameter d generates a least upper bound to fm ■ 



4.4- Charge Renormalization Theory 

Linearized PB theory is valid only for bare macroion charges sufficiently low that 
nonlinear microion screening can be safely neglected. One might conceive of modeling a 
highly charged suspension simply by combining the cell model solution of the nonlinear 
PB equation with the linear-screening approximation for the macroion effective pair 
potential [Eq. (I39l)]. However, without consistently incorporating nonlinear corrections 
to the pair and many-body effective interactions [771 HH] , such a hybrid approximation 
would be uncontrolled. 

A more consistent and computationally practical synthesis of nonlinear screening 
and effective interactions invokes the charge renormalization model of Sec. l2.2l to define 
an effective macroion valence Z within the linear regime. In this scheme, only the N± free 
microions (with — = ZNm), are described by a linear-screening approximation, 
the remaining Z — Z bound counterions (per macroion) renormalizing the bare charge. 
A theory developed previously for salt-free suspensions [^ H5] is here generalized to 
salty suspensions. 

The electrostatic potential around a dressed macroion of effective valence Z and 
association shell thickness 6 is given from Eq. ( 140|1 by 

„k{a+5) p~^^ 

;p{r) = -Z\b --——-^ , r>a + 5, (45) 

1 + K[a + o) r 

with mean value [cf. Eq. fHT]) ] 
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Here h± = N±/[V{1 — f/)] are mean number densities of free microions, fj = 77(1 + S/a)^ 
is the effective volume fraction of the dressed macroions, and k = ^4:TrXB{n+ + n.-) is a 
renormahzed screening constant. From Eqs. ([5]) and f H5l) . the association shell thickness 
is now specified by solving 

ZXb 



C 



(47) 



[1 + k{a + 6)]{a + 6) 

self-consistently for 5 as a function of Z, noting that k depends implicitly on Z and 6. 

The basic distinction between free and bound counterions implies a corresponding 
separation of the volume energy per macroion: 

£v = l^free + /bound • (48) 

The linear-screening approximation is applied only to the free microions, whose grand 
potential (per macroion) is approximated by 

Z^ kXr Z fiA 

T 



i=± 



In 



no 



- 1 



(49) 



1 + k{a + 6) 2 ?2+ + n_ ' 

where x± = N±/Nm- Note that Wfree has the same form as cujj" [Eq. f l30|) ]. except that 
bare parameters are replaced by renormalized parameters and the macroion self energy 
is omitted, being now associated with the bound counterions. 
The bound counterion free energy is approximated by 



bound 



(z-z) 



In 



1 



+ 



2a 



(50) 



the first term on the right side being the ideal-gas free energy of the bound counterions in 
the association shell of volume Vs = {Air / 3) [{a + 6)^ — a^] and the second term accounting 
for the self energy of a dressed macroion, assuming the bound counterions to be localized 
near the macroion surface [r = a). 

Given a bare valence, the effective valence equalizes the chemical potentials of 
counterions in the free and bound phases. This condition is equivalent to minimizing the 
volume energy [Eq. (HSll ] at fixed temperature and mean microion densities [3ll l32l 133] : 



dZ 







(51) 



T,n4 



with Z and 5 related by Eq. fH7j) . The resultant effective valence and corresponding 
shell thickness in turn determine the effective screening constant k. 

The effective pair potential between dressed macroions is given in terms of the 
effective valence and screening constant: 

/ - \ 2 



Z'Xb 



r>2{a + 5) , 



1 + na J r 

from which the macroion free energy can be approximated as [cf. (144|) ] 



mm 

id) 



hs{nm,n±]d) + 2-nnm / dr rVs('^, c?)^eff('^, ^m, 

Jd 



(52) 



(53) 
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The hard-sphere fluid functions, /hs and ^^hs, are calculated from the very accurate 
Carnahan-Starling and Verlet-Weis expressions [88], respectively. In practice, the 
renormalized system parameters (Z, 5, k) must be held constant in the variational 
minimization and in all thermodynamic partial derivatives. 

From the total semigrand potential per macro ion, u = 6^ + fm, the thermodynamic 
pressure finally can be calculated via 



where 



Zinj^ — nJ)k\B 

is the (reduced) pressure due to the entropy and energy of the free microions (bound 
counterions do not contribute directly, since Z and 5 are implicitly fixed in the partial 
derivatives) and 

fiPm = nm + nl^P I 1^ ) (56) 



dn. 



m , 



T,no 

is the contribution from the entropy and effective macroion interactions. The macroion 
pressure also can be obtained from computer simulation using the virial theorem, 
properly generalized to account for the density dependence of the effective pair 
potential [87] . 

Application of charge renormalization theory to bulk suspensions requires setting 
the parameter C in Eqs. (I5|) and (l47l) to define the threshold for renormalization at which 
bare and effective macroion valences begin to differ. Salt-free suspensions afford some 
freedom to set this threshold. In previous applications [HI US], the choice C = 2 yielded 
predictions in close agreement with primitive model simulation data for thermodynamic 
and structural properties of deionized suspensions [89] . For salty suspensions, however, 
the physical requirement that the colon density be non-negative [n-{r) > in Eq. ( l29l) ] 
constrains C < 1. Moreover, the presumed exclusion of colons from the association 
shell (Sec. 12. 2p further dictates that C > 1. Thus, the renormalization threshold is 
uniquely defined by setting C = 1 in the present generalized theory. 



5. Results and Discussion 



We consider the biologically relevant case of aqueous suspensions at room temperature 
{Xb = 0.72 nm) and restrict attention to monovalent microions, justifying neglect 
of microion correlations in the mean-field Poisson-Boltzmann theory. Three different 
implementations of PB theory are here compared: (1) the full nonlinear PB theory 
in a spherical cell model (PB-cell); (2) linearized PB theory in the cell model with 
effective (renormalized) valence (LPB-cell); and (3) linearized PB theory in the effective- 
interaction model (no cell boundaries) with effective valence (LPB-EI). The nonlinear 
PB equation with spherical cell boundary conditions [Eq. ( 12T]) ] was solved numerically 
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by adapting the elegant Mathematica algorithm proposed by Trizac et al. [20] • For the 
linearized PB equation with cell boundary conditions, the electrostatic potential was 
obtained analytically from Eqs. ( I32l) -(l35l) and the pressure numerically from Eq. ( l38l) . 
For the effective-interaction model, the macroion free energy was calculated numerically 
from Eqs. (|52il and (|53ll and the pressure from Eqs. (l5^ -( |56l) . 

Beginning with salt-free suspensions. Figs. [2] and [3] pit PB theory against primitive 
model simulations [7l[8],[9T]. Predictions are compared for the equation of state (pressure 
or osmotic coefficient) at various electrostatic coupling strengths P = Xs/a- Figure [2] 
shows results for Z = 40 and 0.0222 < P < 0.7115, and Fig.[3]for Z = 60 and P = 0.3245. 
All three implementations of PB theory are seen to agree closely with simulation, even 
well above the renormalization threshold [ZT > 5), the parameter regime of highly 
charged latex particles and ionic surfactant micelles. The renormalized linearized 
PB models prove comparable in accuracy to the renormalized jellium model [HI HE]- 
Previous work [^ [87j confirmed the high accuracy of the variational approximation for 
the macroion free energy [Eq. (1531) ] by comparing the LPB-EI model to Monte Carlo 
simulations of the one-component model with the same effective interactions. 

Figure H] compares theoretical predictions for the equation of state of a highly 
deionized suspension at low, but non-zero, salt concentration. For this strongly 
coupled system {ZT ~ 19), the effective macroion valence is significantly lower than 
the bare valence (inset). The same three implementations of PB theory give close 
mutual agreement for the osmotic pressure. In contrast, linearized PB theory with 
bare macroion valence (no renormalization) yields qualitatively different predictions, 
including negative osmotic pressure at lower volume fractions. This unphysical 
behaviour confirms previous findings [lU [45] that charge renormalization can profoundly 
affect thermodynamic properties. When consistently combined with renormalization 
theory, however, the cell and effective-interaction models are comparably accurate in 
describing counterion-dominated suspensions. 

Turning to suspensions in Donnan equilibrium with an electrolyte reservoir, we 
explore the influence of reservoir salt concentration no on osmotic pressure. Figure [5] 
compares predictions for the equation of state of a suspension of highly charged 
macroions with a = 50 nm, Z = 500, over a wide range of concentrations. With 
increasing ionic strength, small relative differences in osmotic pressure for hq < 0.1 mM 
grow to signiflcant deviations for uq > 0.4 mM (Fig. [5^). At concentrations in the 
range 1 mM < rio < 10 mM (Figs. [5b and [5b), a sizeable gap separates the predictions 
of nonlinear PB theory in the cell model and linearized PB theory in the (charge- 
renormalized) effective-interaction model. This gap can be traced to the free energy 
and pressure contributions from microion-induced interactions and correlations between 
macroions, which are naturally included in the effective-interaction model but omitted 
from the cell model [92l [93] . Nevertheless, the salt concentrations predicted by the 
PB-cell and LPB-EI models are in near-exact agreement (inset to Fig. [5t). 

At reservoir concentration no = 1 mM, the cell-model implementation of linearized 
PB theory predicts negative osmotic pressure and bulk modulus over a range of volume 
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fractions. However, since the nonlinear PB-cell model predicts strictly positive pressure 
and bulk modulus [261 EIIj as does the LPB-EI model, this anomalous prediction 
can only be an artifact of the LPB-cell model at intermediate ionic strengths. At 
physiological and higher salt concentrations (no > 100 mM), microion screening so 
weakens effective macroion interactions that accuracy of the cell model is restored for 
the colloidal parameters here studied. 

Finally, we revisit the fundamental issue of thermodynamic phase stability. 
Previous studies based on linearized PB theory [Ml ESI [661 ETJ |68] ^"^^ nonlinear response 
theory [78], without regard for charge renormalization, predicted a spinodal instability 
in deionized suspensions at low (but nonzero) salt concentrations. Separation into 
macroion-rich (liquid) and macroion-poor (vapour) phases was predicted in a parameter 
regime (Z'As/a > 7) that lies slightly beyond the threshold for charge renormalization. 
Subsequent analysis [Ml [25l [261 E3 IM] attributed this unusual instability to an artifact 
of linearization approximations. The theory of Zoetekouw and van Roij [121 133] , which 
incorporates an effective macroion valence, yields a revised prediction of instability 
at much stronger couplings (ZA^/a > 24). This prediction, however, being deduced 
from a mean-field theory, cannot be directly compared with simulations of salt-free 
suspensions [3 [71|9l [121 O [IS] that exhibit macroion aggregation driven by counterion 
correlations. 

In contrast, the present effective-interaction implementation of linearized PB 
theory, surveyed over wide ranges of colloidal parameters (a, Z, rj) and salt 
concentrations, predicts no phase instability for any physically rational renormalization 
threshold (C < 1). Instability emerges only when the theory is pushed beyond 
unphysical renormalization thresholds (C > 1) that imply negative colon densities. The 
absence of thermodynamic instability predicted here is consistent with the renormalized 
jellium model [Ml |35], the cell-model implementation of nonlinear PB theory, which 
admits only a positive bulk modulus [261 IM], and with primitive model simulations for 
weakly correlated microions [95]. The apparent discrepancy between predictions of the 
present approach and that of refs. [121 US] could be related to the different reference 
potentials chosen for linearization of the PB equation. 

6. Conclusions 

Summarizing, a theory of charge renormalization that was previously developed for 
salt-free colloidal suspensions is here generalized to suspensions dialyzed against an 
electrolyte reservoir. The theory incorporates into the one-component model an 
effective macroion valence, defined via a simple thermal criterion. Effective electrostatic 
interactions — one-body volume energy and effective pair potential — are derived from 
Poisson-Boltzmann theory, via a linear-screening approximation, and depend on the 
renormalized parameters. In contrast to the salt-free case, the presence of salt ions 
imposes a unique threshold for charge renormalization. 

When apphed to salt-free, aqueous suspensions with monovalent counterions. 
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the renormalized effective-interaction implementation of linearized PB theory predicts 
osmotic pressures in excellent agreement with the cell-model implementation of 
nonlinear PB theory and with available data from Monte Carlo simulations of the 
primitive model. With increasing salt concentration, predictions of the effective- 
interaction and cell-model implementations increasingly deviate from one another. 
The deviations are attributed to the relatively large pressure contribution originating 
from microion-induced effective interactions and correlations between macroions, which 
are neglected in the cell-model implementation. At physiological and higher salt 
concentrations, however, the macroions are so strongly screened that the relative 
deviations are negligible. We conclude that the cell-model implementation of PB theory, 
although remarkably accurate at relatively low and high salt concentrations, is less 
reliable at intermediate salt concentrations, roughly in the range 0.5 mM< Cg <50 mM. 
More extensive data from controlled experiments and simulations of the primitive model, 
especially for bulk suspensions at significant salt concentrations, are required to test the 
accuracy of competing theories. 

The present approach reveals no indication of a liquid-vapour phase instability. 
Nevertheless, counterion screening in deionized suspensions undoubtedly promotes 
macroion cohesion, as manifested in the density dependence of the volume energy. 
Therefore, effective electrostatic interactions, if not driving, may nevertheless assist an 
instability that is favoured also by other cohesive mechanisms, such as attractive van dcr 
Waals or depletion interactions. The complex interplay between such complementary 
mechanisms for phase separation will be explored in future work. 
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Figure 1. Model of charged colloidal suspension: spherical macroions of radius a and 
point microions dispersed in a dielectric continuum. Strongly associated counterions in 
a spherical shell of thickness 5 renormalize the bare macroion valence Z to an effective 
(reduced) valence Z. 
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Figure 2. Reduced pressure jSp/ntot vs. macroion volume fraction rj, where ntot = 
[Z + (total ion density), of salt-free suspensions with bare macroion valence 

Z = AO and electrostatic coupling strengths (top to bottom) F ~ Xb/o, = 0.0222, 
0.0445, 0.0889, 0.1779, 0.3558, 0.7115. Symbols: Monte Carlo simulations of the 
primitive model [7] (symbol sizes exceed error bars). Curves: Poisson-Boltzmann 
cell model (PB-cell, solid), linearized Poisson-Boltzmann cell model (LPB-cell, dot- 
dashed), and linearized effective-interaction theory (LPB-EI, dashed). For F > 0.1, 
the effective macroion charge is renormalized {Z < Z). 
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Figure 3. Reduced pressure vs. volume fraction 77 of a salt-free suspension with bare 
macroion valence Z — 60 and electrostatic coupling strength F = As/a = 0.3245. 
Symbols: Monte Carlo simulations of the primitive model [8]. Curves: Poisson- 
Boltzmann cell model (PB-cell, solid), linearized Poisson-Boltzmann cell model (LPB- 
cell, dot-dashed), and linearized effective-interaction theory (LPB-EI, dashed). Inset: 
osmotic coefficients over several decades of volume fraction. 
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Figure 4. Reduced osmotic pressure vs. volume fraction of a suspension with 
macroion radius a = 133 nm and bare valence Z = 3500 in Donnan equilibrium with 
a salt reservoir of concentration no = 3.6 /LtM. Predictions of alternative theoretical 
implementations of Poisson-Boltzmann theory are compared: linearized PB theory 
in the effective-interaction (LPB-EI, dashed) and cell-model (LPB-cell, dot-dashed) 
implementations, both incorporating charge renormalization; linearized cell model 
without charge renormalization (LPBO-cell, dotted); nonlinear PB theory in the cell 
model (PB-cell, solid). 
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Figure 5. Reduced osmotic pressure vs. volume fraction of a suspension with 
macroion radius a = 50 nm and bare valence Z = 500 in Donnan equilibrium with 
a salt reservoir of various concentrations: no = 0.1.0.2,0.4 mM (a), 1 mM (b), and 
10 mM (c). Predictions of alternative implementations of Poisson-Boltzmann theory 
are compared: linearized PB theory in the effective-interaction (LPB-EI, dashed) and 
cell-model (LPB-cell, dot-dashed) implementations; nonlinear PB theory in the cell 
model (PB-ccll, solid). Inset to panels (a) and (b): ratio of effective to bare valence 
[predictions of the two LPB models are indistinguishable in panel (b)]. Inset to panel 
(c): ratio of salt concentrations in suspension and reservoir (predictions of the three 
models are indistinguishable on this scale). 



